library(tidymodels)
#> ── Attaching packages ─────────────────────────────────── tidymodels 1.4.1 ──
#> ✔ broom 1.0.9 ✔ recipes 1.3.1
#> ✔ dials 1.4.2 ✔ rsample 1.3.1
#> ✔ dplyr 1.1.4 ✔ tailor 0.1.0
#> ✔ ggplot2 3.5.2 ✔ tidyr 1.3.1
#> ✔ infer 1.0.9 ✔ tune 2.0.0
#> ✔ modeldata 1.5.1 ✔ workflows 1.3.0
#> ✔ parsnip 1.3.3 ✔ workflowsets 1.1.1
#> ✔ purrr 1.1.0 ✔ yardstick 1.3.2
#> ── Conflicts ────────────────────────────────────── tidymodels_conflicts() ──
#> ✖ purrr::discard() masks scales::discard()
#> ✖ dplyr::filter() masks stats::filter()
#> ✖ dplyr::lag() masks stats::lag()
#> ✖ recipes::step() masks stats::step()
data(ames)
ames <- mutate(ames, Sale_Price = log10(Sale_Price))
set.seed(502)
ames_split <- initial_split(ames, prop = 0.80, strata = Sale_Price)
ames_train <- training(ames_split)
ames_test <- testing(ames_split)
lm_model <- linear_reg() %>% set_engine("lm")7 A Model Workflow
在上一章中,我们讨论了parsnip包,该包可用于定义和拟合模型。本章将介绍一个名为模型工作流的新概念。这一概念(以及相应的workflow()对象)的目的是封装建模过程中的主要部分(在1.5节中讨论)。工作流的重要性体现在两个方面:首先,使用工作流概念有助于形成良好的方法体系,因为它是数据分析中估计部分的单一入口;其次,它能让用户更好地组织项目。这两点将在以下各节中详细讨论。
Where Does the Model Begin and End?
到目前为止,当我们使用“模型”这一术语时,指的是将一些预测变量与一个或多个结果变量相关联的结构方程。让我们再以线性回归为例来思考。结果数据表示为yᵢ,其中训练集中有i=1…n个样本。假设模型中使用了p个预测变量xᵢ₁,…,xᵢₚ。线性回归会生成一个模型方程:$ _i = _0 + 1x{i1} + + px{ip} $ 。虽然这是一个线性模型,但它仅在参数上是线性的。预测变量可以是非线性项(例如log(xi))。
关于建模过程的传统思路是,它只包含模型拟合。对于一些简单的数据集,拟合模型本身可能就是整个过程。然而,在拟合模型之前,通常会有多种选择和额外的步骤:
虽然我们的示例模型有p个预测变量,但通常一开始会有多于p个的候选预测变量。通过探索性数据分析或运用领域知识,一些预测变量可能会被排除在分析之外。在其他情况下,可以使用特征选择算法,以数据驱动的方式为模型选择最小的预测变量集。
在某些情况下,重要预测变量的值会缺失。此时,我们不必将该样本从数据集中剔除,而是可以利用数据中的其他值来填补这个缺失值。例如,如果x₁的值缺失,但它与预测变量x₂和x₃相关,那么一种填补方法就可以根据x₂和x₃的值来估计缺失的x₁观测值。
转换预测变量的尺度可能是有益的。如果没有关于新尺度应该是什么的先验信息,我们可以使用统计转换技术、现有数据和一些优化准则来估计合适的尺度。其他转换方法,如主成分分析(PCA),会对一组预测变量进行转换,将它们变成新的特征,用作预测变量。
虽然这些例子都与模型拟合前的步骤有关,但在模型创建后也可能会有一些操作。当创建的分类模型其结果为二元时(例如,event和non-event),通常会使用50%的概率阈值来生成离散的类别预测,这也被称为硬预测。例如,某个分类模型可能会估计某一事件发生的概率为62%。按照通常的默认设置,硬预测结果会是event。然而,该模型可能需要更专注于减少假阳性结果(即,将真正的非事件错误地分类为事件的情况)。实现这一点的一种方法是将阈值从50%提高到某个更高的值。这会提高将新样本判定为事件所需的证据级别。虽然这会降低真阳性率(这是不利的),但它在减少假阳性方面可能会产生更显著的效果。阈值的选择应该利用数据进行优化。这是一个后处理步骤的例子,它对模型的运行效果有着重大影响,尽管它并不包含在模型拟合步骤中。
关注更广泛的建模过程很重要,而不仅仅是拟合用于估计参数的特定模型。这个更广泛的过程包括任何预处理步骤、模型拟合本身以及潜在的后处理活动。在本书中,我们将把这个更全面的概念称为模型工作流(model workflow),并重点介绍如何处理其所有组成部分以生成最终的模型方程。在其他软件中,例如Python或Spark,类似的步骤集合被称为管道。在tidymodels中,“管道(pipeline)”一词已经表示通过管道运算符(例如magrittr包中的%>%或更新的原生|>)链接在一起的一系列操作。为了避免在此语境中使用模糊的术语,我们将与建模相关的一系列计算操作称为工作流。
将数据分析的分析组件整合在一起还有另一个重要原因。后续章节将展示如何准确衡量性能,以及如何优化结构参数(即模型调优)。为了正确量化模型在训练集上的性能,第10章提倡使用重采样方法。要正确做到这一点,分析中任何数据驱动的部分都不应被排除在验证之外。为此,工作流程必须包含所有重要的估计步骤。
举个例子,考虑主成分分析(PCA)信号提取。我们将在第8.4节和第16章中对此进行更详细的讨论;PCA是一种用新的人工特征替代相关预测变量的方法,这些新特征不相关,并且能捕捉原始集合中的大部分信息。这些新特征可以用作预测变量,而最小二乘回归可以用来估计模型参数。
关于模型工作流有两种思考方式。Figure 1 展示了错误的方法:将主成分分析(PCA)预处理步骤视为不属于建模工作流的一部分。
这里的谬误在于,尽管主成分分析(PCA)为生成主成分进行了大量计算,但其运算被假定为不存在任何相关的不确定性。主成分分析的主成分被视为已知的,而且如果不将其纳入模型工作流,就无法充分衡量主成分分析的效果。
Figure 2 展示了一种恰当的方法。
通过这种方式,主成分分析预处理被视为建模过程的一部分。
Workflow Basics
workflows包允许用户将建模对象和预处理对象绑定在一起。让我们再次从Ames数据和一个简单的线性模型开始:
工作流始终需要一个parsnip模型对象:
lm_wflow <-
workflow() %>%
add_model(lm_model)
lm_wflow
#> ══ Workflow ═════════════════════════════════════════════════════════════════
#> Preprocessor: None
#> Model: linear_reg()
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear Regression Model Specification (regression)
#>
#> Computational engine: lm请注意,我们尚未指定此工作流应如何预处理数据:Preprocessor: None。如果我们的模型非常简单,可以使用标准的R公式作为预处理器:
lm_wflow <-
lm_wflow %>%
add_formula(Sale_Price ~ Longitude + Latitude)
lm_wflow
#> ══ Workflow ═════════════════════════════════════════════════════════════════
#> Preprocessor: Formula
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Sale_Price ~ Longitude + Latitude
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear Regression Model Specification (regression)
#>
#> Computational engine: lm工作流有一个fit()方法,可用于创建模型。使用在第6.6节中创建的对象:
lm_fit <- fit(lm_wflow, ames_train)
lm_fit
#> ══ Workflow [trained] ═══════════════════════════════════════════════════════
#> Preprocessor: Formula
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Sale_Price ~ Longitude + Latitude
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#>
#> Call:
#> stats::lm(formula = ..y ~ ., data = data)
#>
#> Coefficients:
#> (Intercept) Longitude Latitude
#> -302.974 -2.075 2.710我们还可以在拟合的工作流上使用predict():
predict(lm_fit, ames_test %>% slice(1:3))
#> # A tibble: 3 × 1
#> .pred
#> <dbl>
#> 1 5.22
#> 2 5.21
#> 3 5.28predict()方法遵循我们在6.3节中为parsnip包描述的所有相同规则和命名约定。模型和预处理器都可以被移除或更新:
lm_fit %>% update_formula(Sale_Price ~ Longitude)
#> ══ Workflow ═════════════════════════════════════════════════════════════════
#> Preprocessor: Formula
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Sale_Price ~ Longitude
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear Regression Model Specification (regression)
#>
#> Computational engine: lm请注意,在这个新对象中,输出显示之前拟合的模型已被移除,因为新公式与之前的模型拟合不一致。
Adding Raw Variables to the workflow()
还有另一种向模型传递数据的接口,即add_variables()函数,它使用类dplyr的语法来选择变量。该函数有两个主要参数:outcomes(结果)和predictors(预测变量)。它们采用类似于tidyverse包的tidyselect后端的选择方法,通过c()来捕获多个选择器。
lm_wflow <-
lm_wflow %>%
remove_formula() %>%
add_variables(outcome = Sale_Price, predictors = c(Longitude, Latitude))
lm_wflow
#> ══ Workflow ═════════════════════════════════════════════════════════════════
#> Preprocessor: Variables
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Outcomes: Sale_Price
#> Predictors: c(Longitude, Latitude)
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear Regression Model Specification (regression)
#>
#> Computational engine: lm也可以使用更通用的选择器来指定预测变量,例如:predictors = c(ends_with("tude"))。一个便利之处在于,任何意外地在predictors参数中指定的结果列都会被悄无声息地移除。这为以下用法提供了便利:predictors = everything()。
当模型拟合时,该声明会将这些数据原封不动地整合到一个数据框中,并将其传递给底层函数:
fit(lm_wflow, ames_train)
#> ══ Workflow [trained] ═══════════════════════════════════════════════════════
#> Preprocessor: Variables
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Outcomes: Sale_Price
#> Predictors: c(Longitude, Latitude)
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#>
#> Call:
#> stats::lm(formula = ..y ~ ., data = data)
#>
#> Coefficients:
#> (Intercept) Longitude Latitude
#> -302.974 -2.075 2.710如果你希望底层建模方法按其通常的方式处理数据,add_variables()会是一个有用的接口。正如我们将在第7.4.1节中看到的,它还能促进更复杂的建模规范。不过,正如我们在下一节中提到的,像glmnet和xgboost这样的模型期望用户从因子预测变量中创建指示变量。在这些情况下,recipe或公式接口通常会是更好的选择。我们将在下一章中介绍更强大的预处理工具——recipe,它也可以添加到工作流程中。
How Does a workflow() Use the Formula?
回顾第3.2节可知,R语言中的公式方法有多种用途(我们将在第8章进一步讨论这一点)。其中之一是将原始数据正确编码为便于分析的格式。这可能包括执行内嵌转换(例如log(x))、创建虚拟变量列、创建交互项或其他列扩展等。然而,许多统计方法需要不同类型的编码:
大多数基于树的模型的R包使用公式接口,但不会将分类预测变量编码为虚拟变量。
某些R包会使用特殊的内联函数来告知模型函数在分析中如何处理预测变量。例如,在生存分析模型中,像
strata(site)这样的公式项会表明列site是一个分层变量。这意味着它不应被当作常规预测变量来处理,并且在模型中没有相应的位置参数估计。有几个R包以base R函数无法解析或执行的方式扩展了公式。在多水平模型(例如混合模型或分层贝叶斯模型)中,诸如
(week | subject)这样的模型项表示列week是一个随机效应,其对于subject列的每个值都有不同的斜率参数估计。
工作流是一种通用接口。当使用add_formula()时,工作流应该如何对数据进行预处理?由于预处理取决于模型,工作流会尽可能模仿底层模型的做法。如果无法做到这一点,公式处理不应对公式中使用的列进行任何操作。让我们更详细地看看这一点。
Special formulas and inline functions
许多多层模型已经采用了lme4包中设计的公式规范。例如,要拟合一个包含受试者随机效应的回归模型,我们会使用下面的公式。其结果是,每个受试者都会有一个针对age的估计截距和斜率参数。
library(lme4)
#> Loading required package: Matrix
#>
#> Attaching package: 'Matrix'
#> The following objects are masked from 'package:tidyr':
#>
#> expand, pack, unpack
library(nlme)
#>
#> Attaching package: 'nlme'
#> The following object is masked from 'package:lme4':
#>
#> lmList
#> The following object is masked from 'package:dplyr':
#>
#> collapse
lmer(distance ~ Sex + (age | Subject), data = Orthodont)
#> Linear mixed model fit by REML ['lmerMod']
#> Formula: distance ~ Sex + (age | Subject)
#> Data: Orthodont
#> REML criterion at convergence: 471.1635
#> Random effects:
#> Groups Name Std.Dev. Corr
#> Subject (Intercept) 7.3912
#> age 0.6943 -0.97
#> Residual 1.3100
#> Number of obs: 108, groups: Subject, 27
#> Fixed Effects:
#> (Intercept) SexFemale
#> 24.517 -2.145标准的R方法无法正确处理这个公式,导致生成一个0行的数据框:
model.matrix(distance ~ Sex + (age | Subject), data = Orthodont)
#> Warning in Ops.ordered(age, Subject): '|' is not meaningful for ordered
#> factors
#> (Intercept) SexFemale age | SubjectTRUE
#> attr(,"assign")
#> [1] 0 1 2
#> attr(,"contrasts")
#> attr(,"contrasts")$Sex
#> [1] "contr.treatment"
#>
#> attr(,"contrasts")$`age | Subject`
#> [1] "contr.treatment"这个特殊公式必须由底层的包代码来处理,而不是标准的model.matrix()方法。即使这个公式可以与model.matrix()一起使用,这仍然会带来问题,因为该公式还指定了模型的统计属性。workflows中的解决方案是使用add_model()补充模型公式——add_variables()明确了预测因子和结果的列名,add_model()中提供实际公式:
library(multilevelmod)
multilevel_spec <- linear_reg() %>% set_engine("lmer")
multilevel_workflow <-
workflow() %>%
# Pass the data along as-is:
add_variables(outcome = distance, predictors = c(Sex, age, Subject)) %>%
add_model(multilevel_spec,
# This formula is given to the model
formula = distance ~ Sex + (age | Subject)
)
multilevel_fit <- fit(multilevel_workflow, data = Orthodont)
multilevel_fit
#> ══ Workflow [trained] ═══════════════════════════════════════════════════════
#> Preprocessor: Variables
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Outcomes: distance
#> Predictors: c(Sex, age, Subject)
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear mixed model fit by REML ['lmerMod']
#> Formula: distance ~ Sex + (age | Subject)
#> Data: data
#> REML criterion at convergence: 471.1635
#> Random effects:
#> Groups Name Std.Dev. Corr
#> Subject (Intercept) 7.3912
#> age 0.6943 -0.97
#> Residual 1.3100
#> Number of obs: 108, groups: Subject, 27
#> Fixed Effects:
#> (Intercept) SexFemale
#> 24.517 -2.145我们甚至可以使用前面提到的生存分析包(survival)中的strata()函数来进行生存分析:
library(censored)
#> Loading required package: survival
parametric_spec <- survival_reg()
parametric_workflow <-
workflow() %>%
add_variables(outcome = c(fustat, futime), predictors = c(age, rx)) %>%
add_model(parametric_spec,
formula = Surv(futime, fustat) ~ age + strata(rx)
)
parametric_fit <- fit(parametric_workflow, data = ovarian)
parametric_fit
#> ══ Workflow [trained] ═══════════════════════════════════════════════════════
#> Preprocessor: Variables
#> Model: survival_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Outcomes: c(fustat, futime)
#> Predictors: c(age, rx)
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Call:
#> survival::survreg(formula = Surv(futime, fustat) ~ age + strata(rx),
#> data = data, model = TRUE)
#>
#> Coefficients:
#> (Intercept) age
#> 12.8734120 -0.1033569
#>
#> Scale:
#> rx=1 rx=2
#> 0.7695509 0.4703602
#>
#> Loglik(model)= -89.4 Loglik(intercept only)= -97.1
#> Chisq= 15.36 on 1 degrees of freedom, p= 8.88e-05
#> n= 26Creating Multiple Workflows at Once
在某些情况下,数据需要多次尝试才能找到合适的模型。例如:
对于预测模型,建议评估多种不同的模型类型。这需要用户创建多个模型规格。
模型的序贯检验通常从一组扩展的预测变量开始。这种“全模型”会与一系列相同的模型进行比较,这些模型依次移除了每个预测变量。通过基本的假设检验方法或经验验证,可以分离并评估每个预测变量的影响。
在这些情况下以及其他一些情况中,根据不同的预处理程序集和/或模型规格创建大量工作流可能会变得繁琐或繁重。为了解决这个问题,workflowset包会创建工作流组件的组合。一系列预处理程序(例如公式、dplyr选择器,或下一章将讨论的特征工程recipe对象)可以与一系列模型规格相结合,从而生成一套工作流。
举个例子,假设我们想重点关注Ames数据中房屋位置的不同表示方式。我们可以创建一组公式来捕捉这些预测变量:
location <- list(
longitude = Sale_Price ~ Longitude,
latitude = Sale_Price ~ Latitude,
coords = Sale_Price ~ Longitude + Latitude,
neighborhood = Sale_Price ~ Neighborhood
)可以使用workflow_set()函数将这些表示形式与一个或多个模型进行组合。我们将仅使用之前的线性模型规格来进行演示:
library(workflowsets)
location_models <- workflow_set(preproc = location, models = list(lm = lm_model))
location_models
#> # A workflow set/tibble: 4 × 4
#> wflow_id info option result
#> <chr> <list> <list> <list>
#> 1 longitude_lm <tibble [1 × 4]> <opts[0]> <list [0]>
#> 2 latitude_lm <tibble [1 × 4]> <opts[0]> <list [0]>
#> 3 coords_lm <tibble [1 × 4]> <opts[0]> <list [0]>
#> 4 neighborhood_lm <tibble [1 × 4]> <opts[0]> <list [0]>
location_models$info[[1]]
#> # A tibble: 1 × 4
#> workflow preproc model comment
#> <list> <chr> <chr> <chr>
#> 1 <workflow> formula linear_reg ""
extract_workflow(location_models, id = "coords_lm")
#> ══ Workflow ═════════════════════════════════════════════════════════════════
#> Preprocessor: Formula
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Sale_Price ~ Longitude + Latitude
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#> Linear Regression Model Specification (regression)
#>
#> Computational engine: lm工作流集主要设计用于处理重采样,这在第10章中会进行讨论。option列和result列必须填充从重采样中得到的特定类型的对象。我们将在第11章和第15章中更详细地演示这一点。
与此同时,让我们为每个公式创建模型拟合,并将它们保存在一个名为fit的新列中。我们将使用基本的dplyr和purrr操作:
location_models <-
location_models %>%
mutate(fit = map(info, ~ fit(.x$workflow[[1]], ames_train)))
location_models
#> # A workflow set/tibble: 4 × 5
#> wflow_id info option result fit
#> <chr> <list> <list> <list> <list>
#> 1 longitude_lm <tibble [1 × 4]> <opts[0]> <list [0]> <workflow>
#> 2 latitude_lm <tibble [1 × 4]> <opts[0]> <list [0]> <workflow>
#> 3 coords_lm <tibble [1 × 4]> <opts[0]> <list [0]> <workflow>
#> 4 neighborhood_lm <tibble [1 × 4]> <opts[0]> <list [0]> <workflow>
location_models$fit[[1]]
#> ══ Workflow [trained] ═══════════════════════════════════════════════════════
#> Preprocessor: Formula
#> Model: linear_reg()
#>
#> ── Preprocessor ─────────────────────────────────────────────────────────────
#> Sale_Price ~ Longitude
#>
#> ── Model ────────────────────────────────────────────────────────────────────
#>
#> Call:
#> stats::lm(formula = ..y ~ ., data = data)
#>
#> Coefficients:
#> (Intercept) Longitude
#> -184.396 -2.025我们在这里使用了一个purrr函数来遍历我们的模型,但有一种更简单、更好的方法来拟合工作流集合,这将在第11.1节中介绍。
总的来说,工作流集还有很多内容值得探讨!虽然我们在这里介绍了基础知识,但工作流集的细微差别和优势要到第15章才会详细阐述。
Evaluating the Test Set
假设我们已经完成了模型开发,并确定了最终模型。有一个便捷函数叫做last_fit(),它会用整个训练集对模型进行拟合,并用测试集对其进行评估。
以lm_wflow为例,我们可以将模型以及初始的训练/测试集拆分传递给该函数:
final_lm_res <- last_fit(lm_wflow, ames_split)
final_lm_res
#> # Resampling results
#> # Manual resampling
#> # A tibble: 1 × 6
#> splits id .metrics .notes .predictions
#> <list> <chr> <list> <list> <list>
#> 1 <split [2342/588]> train/test split <tibble [2 × 4]> <tibble> <tibble>
#> # ℹ 1 more variable: .workflow <list>请注意,last_fit()接受的数据是已分割的,而不是原始数据框。该函数使用此分割好的数据集来生成用于最终拟合和评估的训练集和测试集。
.workflow列包含拟合的工作流,可以使用以下代码从结果中提取出来:
fitted_lm_wflow <- extract_workflow(final_lm_res)同样地,collect_metrics()和collect_predictions()分别提供了获取性能指标和预测结果的途径。
collect_metrics(final_lm_res)
#> # A tibble: 2 × 4
#> .metric .estimator .estimate .config
#> <chr> <chr> <dbl> <chr>
#> 1 rmse standard 0.164 pre0_mod0_post0
#> 2 rsq standard 0.189 pre0_mod0_post0
collect_predictions(final_lm_res) %>% slice(1:5)
#> # A tibble: 5 × 5
#> .pred id Sale_Price .row .config
#> <dbl> <chr> <dbl> <int> <chr>
#> 1 5.22 train/test split 5.02 2 pre0_mod0_post0
#> 2 5.21 train/test split 5.39 4 pre0_mod0_post0
#> 3 5.28 train/test split 5.28 5 pre0_mod0_post0
#> 4 5.27 train/test split 5.28 8 pre0_mod0_post0
#> 5 5.28 train/test split 5.28 10 pre0_mod0_post0我们将在16.6节中进一步了解last_fit()的实际应用以及如何再次使用它。
在使用验证集时,last_fit()有一个名为add_validation_set的参数,用于指定最终模型是仅在训练集上训练(默认设置),还是在训练集和验证集的组合上训练。
Chapter Summary
在本章中,你了解到建模过程不仅仅包括估计将预测变量与结果联系起来的算法参数。这个过程还包括预处理步骤以及模型拟合后的操作。我们引入了一个名为模型工作流的概念,它可以捕捉建模过程中的重要组成部分。多个工作流也可以在一个工作流集合中创建。last_fit()函数便于将最终模型拟合到训练集并使用测试集进行评估。
对于Ames数据集,其使用的相关代码如下:
library(tidymodels)
data(ames)
ames <- mutate(ames, Sale_Price = log10(Sale_Price))
set.seed(502)
ames_split <- initial_split(ames, prop = 0.80, strata = Sale_Price)
ames_train <- training(ames_split)
ames_test <- testing(ames_split)
lm_model <- linear_reg() %>% set_engine("lm")
lm_wflow <-
workflow() %>%
add_model(lm_model) %>%
add_variables(outcome = Sale_Price, predictors = c(Longitude, Latitude))
lm_fit <- fit(lm_wflow, ames_train)